Song et al. BMC Genomics 2014, 15:297 
http://www.bionnedcentral.conn/1471 -21 64/1 5/297 



RESEARCH ARTICLE Open Access 



The phenotypic predisposition of the parent in F1 
hybrid is correlated with transcriptome 
preference of the positive general combining 
ability parent 

Gaoyuan Song, Zhibin Guo, Zhenwei Liu, Xuefeng Qu, Daiming Jiang, Wei Wang, Yingguo Zhu 
and Daichang Yang"" 



Abstract 

Background: Sprague and latum (1942) introduced the concepts of general combining ability (GCA) and specific 
combining ability (SCA) to evaluate the breeding parents and Fl hybrid performance, respectively. Since then, the 
GCA was widely used in cross breeding for elite parent selection. However, the molecular basis of GCA remains to 
unknown. 

Results: We studied the transcriptomes of three varieties and three Fl hybrids using RNA-Sequencing. Transcriptome 
sequence analysis revealed that the transcriptome profiles of the Fls were similar to the positive GCA-effect parent. 
Moreover, the expression levels of most differentially expressed genes (DEGs) were equal to the parent with a positive 
GCA effect. Analysis of the gene expression patterns of gibberellic acid (GA) and flowering time pathways that 
determine plant height and flowering time in rice validated the preferential transcriptome expression of the parents 
with positive GCA effect. Furthermore, H3K36me3 modification bias in the Pseudo-Response Regulators (PRR) gene family 
was observed in the positive GCA effect parents and demonstrated that the phenotype and transcriptome bias in the 
positive GCA effect parents have been epigenetically regulated by either global modification or specific signaling 
pathways in rice. 

Conclusions: The results revealed that the transcriptome profiles and DEGs in the Fls were highly related to 
phenotype bias to the positive GCA-effect parent. The transcriptome bias toward high GCA parents in Fl hybrids 
attributed to H3K36me3 modification both on global modification level and specific signaling pathways. Our results 
indicated the transcriptome profile and epigenetic modification level bias to high GCA parents could be the molecular 
basis of GCA. 
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Background 

Selecting elite parents is of paramount importance in 
cross- and hybrid-breeding programs. The selection of 
parents from a phalanx of inbred lines, however, is ex- 
tremely laborious and time-consuming and can be ran- 
dom. Adding to this complexity, parents with excellent 
agronomic traits do not always pass those traits on to 
their progeny. To evaluate breeding parents, Sprague 
and Tatum (1942) introduced the concepts of general 
combining ability (GCA) and specific combining ability 
(SCA), which allow the study and comparison of the 
performances of inbred lines in hybrid combination. 
GCA is used to designate the average performance of an 
inbred line in hybrid combination, and SCA is used to 
designate those cases in which certain combinations do 
relatively better or worse than would be expected on the 
basis of the average performance of the lines involved 
[1]. GCA and SCA effects have been successfully used as 
criteria to evaluate elite parents in conventional cross- 
breeding and the performance of hybrid combinations. 

Since the introduction of the GCA concept in 1942, 
very limited genetic studies on GCA have been con- 
ducted even though GCA is widely used in breeding 
programs for evaluation of the parents in early genera- 
tions [2-5]. A few quantitative genetic loci with GCA 
have been identified recently. Qu et al analyzed the 
QTLs of 10 agronomic traits for GCA using recombin- 
ant inbred line (RIL) populations with three testers in 
three testcross populations and a backcross recombinant 
inbred line (BCRIL) population of rice [6] . They detected 
a large number of additive effects of QTL^^^ loci. Qi 
et al found that several genetic loci responding for 
GCA and SCA for five yield-related traits using a set of 
testcrosses with introgression lines (ILs) of maize under 
different environmental conditions. Total of 56 signifi- 
cant QTL^^^ loci have been mapped [7]. These studies 
have revealed that GCA effects, like the traits, are genet- 
ically controlled. At the molecular level, however, how 
the phenotypes associated with GCA effects are passed 
on to the Fl hybrids remains unknown. 

High-throughput genome-wide analysis approaches 
such as microarray analysis and next generation sequen- 
cing have been used in the study of phenotypes such as 
hybrid vigor [8-10]. Genome-wide gene expression pro- 
files related to heterosis have been studied [11-14]. Stupar 
et al. studied the gene expression profiles between 
maize Fl hybrids and their parents, and approximately 
75% of the differentially expressed genes showed addi- 
tive expression patterns between Fl hybrids and par- 
ents [13]. They studied the genetic diversity and 
transcriptional variation with different maize hybrids 
and found that the genetic diversity was correlated with 
transcriptional variation, and little (less than 1%) of the 
gene expression in Fl hybrids was outside the parental 



range [13]. Wei et al studied the gene expression profile 
between the super hybrid rice LYP9 and its parents and 
suggested that the differentially expressed genes might 
correlate with heterosis. Moreover, Riedelsheimer et al 
have analyzed metabolic pathways in maize hybrids and 
were able to predict GCA scores using genome-wide asso- 
ciation analysis [15]. They have used 285 crosses derived 
from the diverse inbred lines of maize with two testers 
and predicted their combining ability for seven biomass- 
and bioenergy related traits using 56,110 single nucleotide 
polymorphisms (SNPs) and 130 metabolites. The predic- 
tion accuracies were from 0.72 to 0.81 for SNPs and 0.60 
to 0.80 for metabolites. A genome-wide analysis of the as- 
sociation between GCA and gene expression profiles has 
not been investigated. The contribution of the elite par- 
ents with GCA effects to their offspring also remains 
unexplored. 

In the present study, we constructed a diallele crossing 
population with five rice parents and 10 Fl hybrids. We 
evaluated the GCA effects of three agronomic traits re- 
lated to grain yield, heading date, plant height and grain 
number from five elite rice varieties. The results showed 
that 93-11 and Teqing (TQ) have positive GCA effects, 
and Guangluai 4# (GL), Aijiaonante (AJ) and Zhenshan 
97 (ZS) showed negative GCA effects. The analysis of 
the transcriptome profiles of the leaves from three Fl 
hybrids, GL x 93-11, GL x TQ and 93-11 x TQ, and their 
parents revealed that transcriptome profiles were corre- 
lated to the positive-GCA-effects parent, showing obvi- 
ous parental bias. Further analysis found that the 
expression levels of most of the DEGs were obviously 
biased towards to the positive GCA effect parent. Ana- 
lysis of the gene expression patterns of gibberellic acid 
(GA) and flowering-time signaling pathways for plant 
height and flowering time validated the transcriptome 
bias to the positive GCA effect parents. Our results also 
indicated that the H3K36me3 modification showed bias 
to the positive GCA effect parent and demonstrated that 
phenotypes of the positive GCA effect parents were at- 
tributed to transcriptome bias. These results will be 
helpful in understanding the parental gene contribution 
to Fl hybrids and the molecular basis of GCA. 

Results 

Evaluation of general combining ability of elite rice 
varieties 

Five inbred rice varieties representing different breeding 
objectives from the 1970s to the present were chosen. 
Ten Fl hybrids were obtained from five varieties using a 
diallele crossing design (Additional file 1: Figure SI and 
Additional file 2: Figure S2). Three important agronomic 
traits related to grain yield (that is, plant height, heading 
date and grains per panicle) were evaluated based on the 
phenotypes of Fl hybrids from the diallele crossing 
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populations (Figure lA, B and C). The results revealed 
that the GCA scores for the two rice varieties 93-11 and 
TQ were positive, whereas negative GCA scores were 
observed for the varieties GL, ZS and AJ (Figure ID). 
Further analysis revealed that the phenotypes of the 
three agronomic traits in the Fl hybrids were always 
similar to or higher than 93-11 and TQ, when either 
was used as one parent in the cross (Figure lA-C and 
Additional file 1: Figure SI). In contrast, the phenotypes 
of Fl hybrids derived from GL (Figure 1), AJ and ZS 
(Additional file 1: Figure SI and Additional file 2: Figure 
S2) were significantly different compared to 93-11 and 
TQ. The parental phenotypes were significantly corre- 
lated with the Fl hybrids crossed to 93-11 and TQ (r > 
0.97, P < 3.80E-07) but less significantly with the Fl hy- 
brids resulting from crosses to GL, ZS and AJ (r = 0.81, 
P < 1.33E-03) (Table 1). These results indicated that the 
GCAs for the three agronomic traits relevant to grain 
yield in parents 93-11 and TQ were significantly higher 
than those of the varieties GL, AJ and ZS. 

The phenotype bias toward the positive GCA effect 
parent is attributable to transcriptome bias 

As shown in Figure 1 and Additional file 3: Figure S3, 
GL had an earlier heading date (68.4 days on average) 
and had a shorter plant height (88.9 cm on average), 
whereas TQ and 93-11 had later heading dates (90.1 
and 93.7 days on average) and had a taller plant height 
(113.3 and 117.2 cm on average, respectively). The heading 



Table 1 The correlation of phenotypes between each 
parent and their F1s 



Parents 


r 


P value 


ZS 


0.81 


1.33E-03 


AJ 


0.81 


1.27E-03 


GL 


0.81 


1.26E-03 


93-11 


0.97 


3.80E-07 


TQ 


0.99 


2.06E-10 



dates (90.5 and 109.5 days on average, respectively) and 
plant heights (113.9 and 117.0 cm on average, respect- 
ively) in the Fl hybrids derived from TQ and 93-11 
were obviously biased toward TQ or 93-11. To under- 
stand the relationship between the Fl phenotype and 
the positive GCA effect parent, we conducted a tran- 
scriptome profile analysis of the three parents, GL, 93-11 
and TQ, and their three Fl hybrids, GLx 93-11, GL x 
TQ and 93-11 x TQ, using RNA sequencing technology 
(Additional file 3: Figure S3) [11]. A cluster analysis in- 
dicated that the transcriptomes of GL x TQ and GL x 
93-11 are similar to those of the positive GCA effect 
parents, TQ and 93-11, respectively, and are signifi- 
cantly different from those of the negative GCA effect 
parent GL (Figure 2A and 2B). The transcriptome pro- 
file of 93-11 X TQ was more similar to 93-11 that had 
the higher positive GCA effect than to TQ (Figure 2C). 
The transcriptome similarity of the Fl hybrids to either 
the 93-11 or TQ parent was consistent with the 
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Figure 2 Cluster analysis of global transcriptome of the three F1 hybrids and their parents. Clusters of transcriptomes between GL, GLx 
TQ and TQ (A); GL, GLx 93-11 and 93-11 (B) and 93-11, 93-11 xTQ and TQ (C). The percentage of DEGs between Fl hybrids and the parents, 
GLxTQ (D), GLx 93-11 (E) and 93-11 xTQ (F). 



pheno types of the three traits in 93-11 and TQ. To ex- 
plore the transcriptome profiles of the parents in the Fl 
hybrids, we analyzed the gene expression level in Fl hy- 
brids. The results showed that 76.0% of the gene expres- 
sion profiles in GL x TQ were similar to GL and 83.7% 
were similar to TQ (Figure 2D). Analogous results were 
found for the GL x 93-11 and 93-11 x TQ (i.e., 84% of 
the gene expression profiles in GL x 93-11 were similar 
to those in 93-11 and 67.1% of genes bias to GL) (Figure 2E). 
Of 90.2% expressed genes in 93-11 x TQ were similar 
to those in 93-11 and 82.6% of the genes bias to TQ 
(Figure 2F). Therefore, the phenotypes of the Fls de- 
rived from the positive GCA effect parents are corre- 
lated with a transcriptome bias in the Fls. 

Differentially expressed genes in the Fl hybrids are those 
preferentially expressed in the positive GCA effect parent 

Differentially expressed genes (DEGs) have been recog- 
nized to play important biological functions in heterosis 
[16-18]. To explore whether DEGs contribute to the 
positive GCA effect parent in the Fl hybrids, we ana- 
lyzed the DEG profiles from three Fl hybrids. A total of 
34,486-35,718 genes were expressed in the Fl hybrids 
and their parents (Figure 3A). We found that 22.1% to 
37.0% genes were differentially expressed in the Fl 



hybrids comparing GL x TQ, GL x 93-11 and 93-11 x TQ 
with their parents, respectively (fold change >2.0, P < 
0.05) Moreover, 75.5% to 84.1% of the DEGs between the 
Fl hybrids and the parents had similar expression levels 
to one parent (Figure 3B). Of these, 63.4%, 79.5% and 
69.6% of the gene expression levels were similar to the 
positive GCA effect parents TQ and 93-11 in GL x TQ, 
GLx 93-11 and 93-11 x TQ, respectively (Figure 3C). 
The expression levels of the remaining genes were simi- 
lar to the negative or lower GCA effect parents GL and 
TQ (Figure 3C). 7.9% to 19.7% of the DEGs expression 
level in the Fl hybrids were out of the parental ranges, 
whereas 4.8% to 10.2% of the DEGs showed mid-parent 
expression levels, but showing significant difference to 
both parents (Figure 3B). These results reveal that the 
expression levels of the majority of the DEGs were simi- 
lar to those of the positive GCA parent. 

Expression patterns of the genes in the regulatory 
pathways controlling flowering time and plant height in 
rice 

To further support the observed correlation between the 
positive GCA effect pheno type in the parent and tran- 
scriptome bias in the Fl, we chose the genes in two well- 
known signaling pathways controlling flowering time 
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A. 





Total 


Non-DEGs 


% 


DEGs 


% 


GLXTQ vs Parents 


35,097 


23.742 


67.6% 


11,355 


32.4% 


GLX93-1 1 vs Parents 


35,718 


22,493 


63.0% 


13,225 


37.0% 


93-1 1XTQ vs Parents 


34,486 


26,870 


77.9% 


7,616 


22.1% 




GLXTQ 



GLX93-1 1 



93-1 1XTQ 




Figure 3 Global differentially expressed genes between the F1 hybrids and their parents. (A) The proportions of DEGs and non-DEGS were 
detected in tliree Fl liybrids; (B) categorization of different types of DEGs between Fl liybrids and parents. Equal to parent means tlie expression 
level difference in the Fl hybrids exists either one parent; over parents means the expression level in Fl the hybrids is out of the ranges of both 
parents; between parents means the difference of the expression level in Fl hybrids exists both parents, but not out of the parent range; and 
(C) Parent bias of the DEGs that were expressed at a level equal to one parent. 



(Additional file 4: Table SI) and plant height in rice. The 
expression oi EARLY HEADING DATE 1 (Ehdl) smd EARLY 
HEADING DATE2 {Ehd2) and repression of HEADING 
DATEl (Hdl), HEADING DATESa (HdSa) and RICE 
FLOWERING LOCUS Tl (RFTl) are activated in a cas- 
cade to initiate rice flowering [19-25]. Expression of the 
Rice Pseudo-Response Regulators 1 {OsPRRl) and Grain 
number, plant height and heading date? {GhdT) involves 
a delay in the flowering time and an increase in the 
plant height and grain number [26-28]. We found that 
the expression levels of OsPRRl, Hdl and Ghd7 were 
very high in TQ, 93-11, GL x TQ, GL x 93-11 and 93- 
11 X TQ versus GL (Table 2). By contrast, the expression 
levels of RFTl, HdSa, Ehdl and Ehd2 were higher in GL 
(Table 2). These results showed that the expression pat- 
terns of the flowering regulation genes were consistent 
with transcriptome bias towards the positive GCA effect 
parents in the Fl hybrids. 

The expression patterns of the other set of genes in- 
volved in gibberellic acid (GA) metabolism and the 
signaling pathway that controls plant height in rice 
were also analyzed between Fls and parents [29-32] 
(Additional file 5: Table S2). In plants, bioactive GAs are 
synthesized from the precursor geranylgeranyl diphos- 
phate by ent-copalyl diphosphate synthase (CPS) and ent- 
kaurene synthase (KS), followed by ent-kaurene oxidase 
(KO) and ent-kaurenoic acid oxidase (KAO). At the final 
stage, GA2o-oxidase (GA20ox) and GAs-oxidase (GA3ox) 
catalyze the conversion of GA53/GA12 and GA1/GA4. The 
bioactive GAs and the precursors elongated the upper- 
most inter node and were deactivated by GA2- oxidase 
(GA2ox) and (EUI) [30]. As shown in Table 2, three genes, 
ent-copalyl diphosphate synthase {OsCPSl), ent-kaurenoic 



acid oxidase (OsKAO) and GA2-oxidase 6 {OsGA2ox6) 
had the higher expression levels in TQ, 93-11, GL x TQ, 
GL X 93-11 than in GL. In contrast, expression levels of 
the other three genes, GA3ox2 and GA20ox2 were less in 
TQ, 93-11, GL X TQ, GL x 93-11 and 93-11 x TQ than 
in GL. The expression profiles of those genes exactly 
matched the feedback and feed-forward regulation 
mechanism of GAs synthesis [31,32]. Our results are 
consistence with the results reported in wheat hybrids 
[33,34]. The results again demonstrated that the expres- 
sion patterns of the genes corresponding to plant height 



Table 2 Expression level of flowering- and GA 
metabolism-related genes 



Gene name 




Expression level (RPKM) 








GL 


TQ 


93-11 


GLXTQ 


GLX93-1 1 


93-1 1XTQ 


Flowering pathway 












OsPRRl 


35.49 


92.53 


122.96 


122.26 


89.91 


138.95 


RFTl 


2.28 


0.16 


0.43 


1.45 


0.48 


0.20 


Hd3o 


1.33 


0.00 


0.00 


0.21 


0.00 


0.00 


Hdl 


0.57 


3.27 


1.30 


2.17 


5.82 


1.62 


GHd7 


0 


0 


0.30 


0 


0.64 


0.18 


EHd2 


2.05 


0.83 


0.22 


1.08 


0.81 


0.32 


EHdl 


2.63 


0 


0 


0 


0 


0 


GA metabolism pathway 










OsCPSl 


0.23 


0.46 


0.71 


0.92 


0.50 


0.67 


OsKAO 


0.07 


0.61 


0.66 


0.95 


1.17 


0.61 


OsGA3ox2 


2.12 


0.10 


0.22 


0.80 


0.46 


0.92 


OsGA20ox2 


12.92 


3.15 


4.84 


9.36 


5.78 


5.51 


OsGA2ox6 


0.82 


5.43 


1.12 


4.54 


8.37 


7.95 
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exhibited a bias toward the positive GCA effect parents. 
Taken together, the data validated the transcriptome 
bias toward the positive GCA effect parents through in- 
dividual metabolism pathways for plant height and 
heading date in rice and further demonstrated that the 
phenotype bias to the positive GCA effect parent in Fl 
hybrids is due to the transcriptome bias toward the 
positive GCA effect parents. 

Transcriptome bias toward the positive GCA effect parent 
regulated by H3K36me3 modifications 

Previous studies have described the histone modifications 
involved in the regulation of the transcriptome [35-38]. 
Previous studies have indicated that the trimethylation of 
histone H3 on lysine 27 (H3K27me3) in gene body re- 
presses gene expression [10,39,40] and trimethylated his- 
tone H3 on lysine 36 (H3K36me3) in the gene body 
activated gene expression [41-43]. To further explore the 
mechanism of transcriptome bias to the positive GCA ef- 
fect parent in Fl hybrids, we analyzed the patterns of 
methylation at H3K36me3 and H3K27me3 in Fl hybrids. 
Genome-wide histone modifications of H3K36me3 and 
H3K27me3 in GLx 93-11, GL x TQ and their parents, 
GL, 93-11 and TQ, were analyzed with the aid of a DNA 
library prepared after chromatin immunoprecipitation 
(Chip). The results showed 18.2%-18.3% genes were over- 
lapped between H3K36me3 and H3K27me3 modifica- 
tions. However, an obvious parental bias of H3K36me3 
modification was detected. In GLx 93-11, 82.0% of the 




GLX93-11 vs 



genes with H3K36me3 modifications were similar to 
93-11 (FDR > 0.001), whereas 73.4% were similar to GL 
(FDR > 0.001) (Figure 4A). 93.0% of H3K36me3 modifi- 
cations in GL X TQ were found to prefer TQ (FDR > 
0.001), and 75.0% of H3K36me3 preferred GL (FDR > 
0.001). For H3K27me3 modifications, 69.8% of the 
genes in GL x 93-11 with modifications were similar to 
GL (FDR > 0.001), and 65.2% were biased to 93-11 (FDR > 
0.001) (Figure 4B). In GL x TQ, 90.7% genes with the 
H3K27me3 modifications preferred TQ (FDR > 0.001), 
whereas only 62.5% of H3K27me3 modifications were 
biased to GL (FDR > 0.001) (Figure 4B). There was no 
significant correlation of H3K27me3 between positive 
GCA effect and negative GCA effect parents in two Fl 
hybrids. Taken together, the results indicated that 
H3K36me3 modifications exhibited a bias toward the 
positive GCA effect parent in the Fl hybrids, which sug- 
gested that epigenetic regulation could involve regula- 
tion of the phenotypes and transcriptome bias toward 
positive GCA effect parents. 

Furthermore, we analyzed the modification patterns of 
Cryptochrome 2 {CRY2), Pseudo-Response Regulator 37 
and 95 {PRR37 and PRR9S) in Fl hybrids (Figure 5). 
These genes are involved in the regulation of rice flower- 
ing time and adaptability. CRY2 is a photolyase-like 
blue-light receptor that mediates light responses in 
plants via interaction with the CIBl (cryptochrome- 
interacting basic-helix-loop-helix) protein to promote 
Ci?r2-dependent floral initiation [44]. The PRR37 gene 



■ DMG with GL 

■DMG with 93-11 

DMG with both 
parents 

■ Non DMG 




' DMG with GL 



' DMG with TQ CO 
0 

DMG with both ^ 
parents CO 

I Non DMG ^0 
X 
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GLXTQ vs parents 




GLX93-1 1 vs parents GLXTQ vs parents 

Figure 4 The epigenetic modification level in Fl hybrids compared to each parent. (A) Percentages of genes with H3K36nne3 modification 
levels similar to those in each parent in Fl hybrids. (B) Percentages of genes with H3K27me3 modification levels similar to those of each parent 
in Fl hybrids. 
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Figure 5 H3K36me3 modification level of some flowering pathway genes between F1 hybrids and parents. (A) H3K36me3 modification 
level between GLx93-l 1 and the parents GL, 93-1 1. (B) H3K36me3 modification level between GLxTQ and the parents GL, TQ. 



is involved in the down-regulation of Hd3a gene expres- 
sion to suppress flowering under long-day conditions. 
The varieties harboring nonfunctional alleles of PRR37 
flower extremely early under natural, long-day condi- 
tions [45]. We found that the H3K36me3 modification 
levels of three genes in GLx 93-11 or GL x TQ biased 
toward the positive GCA effect parents, TQ or 93-11 
(Figure 5). The results are consistent with the heading 
date phenotype bias toward high GCA effect parents in 
the Fl. Again, our data demonstrated that the phenotype 
and transcriptome biased toward the positive GCA effect 
parents were epigenetically regulated by both global 
modification and specific signaling pathways in rice. 

Discussion 

The GCA concept has been studied for more than 
70 years [1]. It has being successfully applied in crop 
and livestock breeding for the evaluation of parent per- 
formance [3-5]. However, in spite of this, the genetic 
analysis of the GCA is very limited [6,7], and the mo- 
lecular mechanism of the GCA has not been docu- 
mented. In the present study, we found that certain 
phenotypes in rice Fl hybrids were always biased toward 
the positive GCA effect parents in the traits tested (r > 
0.97, P <3.80E-07). Further studies showed that the 
transcriptome profiles in Fl hybrids were the same as 
phenotypes that are biased to the parents with positive 
GCA effect. The expression patterns of individual genes 
corresponding to plant height and flowering time dem- 
onstrated that the phenotype bias toward the positive 
GCA effect parent is attributed to the transcriptome 
profile and the specific epigenetic modification. Our data 
revealed that the phenotype bias in Fl hybrids could be 
attributed to transcriptome and H3K36me3 modification 
bias. Our findings provide molecular clues as to why and 
how the positive GCA effect parent has been widely 
used for crop and livestock improvement in breeding 
programs for half a century. 

Although GCA is widely recognized and applied to crop 
and livestock improvement by breeders, it is based largely 
on breeding experiences, and therefore, it highly is unpre- 
dictable. Our findings provide molecular evidence of GCA 



through transcriptome and epi-genome analysis. This 
study first revealed that the transcriptome and epigenome 
in rice Fl hybrids are similar to the parent with the posi- 
tive GCA performance, but significantly different from the 
negative GCA effect parents. We also found that the posi- 
tive GCA effect parents showed better performance with 
regard to agronomic traits in Fl hybrids compared to the 
negative GCA effect parents. In the previous studies, a 
large number of genetic loci with dominance effect were 
detected in Fl hybrids [46-48]. eQTL analysis revealed 
that the gene expression regulation was an complicated 
regulation networks [49-51]. So, the phenomena of the 
phenotypes and the transcriptomes biased to high GCA 
parents in Fl hybrids might be consequence of the accu- 
mulation of the loci with dominance effect in the elite par- 
ents through artificial selection. Our results suggested that 
the transcriptomes profile in Fl hybrids could be opti- 
mized during elite parent selection, which give rise to elite 
agronomic traits, such as increase of biomass and grain 
yield. During the parent improvement process, these al- 
leles corresponding to the agronomic traits satisfying dif- 
ferent breeding objectives had been selected or kept. For 
instance, a representative elite variety of GL that met the 
breeding objectives of more panicle and less grain number 
per panicle had been developed in the 1970s in China. 
Later, an ideal plant type that has less panicle numbers 
and more spikelet per panicle was raised, and the breeding 
objectives of higher grain yields required high biomass 
that was due to later heading date and higher plant height. 
To meet the requirements of the breeding objective, the 
alleles corresponding to late heading date and higher plant 
height were accumulated and optimized through artificial 
selection pressure. Therefore, the modern elite parent 
93-11 has a longer flowering time and taller plant height 
that produces higher biomass resulting in higher grain 
yield. Consequently, the alleles corresponding to later flow- 
ering time and taller plant height were selected and main- 
tained in the genome of the elite parent. Our results again 
demonstrated the signiflcance of selecting elite parents with 
excellent GCA performance in crossbreeding program. 

Several studies have indicated that the mechanism of 
the aUelic-specific expression could be regulated by 
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epigenetic modification, including DNA methylation and 
various histone modifications [52,53]. Previous studies 
have detected maternal alleles for some imprinted genes 
were hypomethylated in the endosperm, which con- 
tributed to the parent-of-origin expression pattern in 
reciprocal hybrids in maize [52,53]. However, very few 
parent-of-origin effects were detected in plant vegetative 
tissues [54]. Our global gene expression analysis has not 
found any parent-of-origin genes, but almost of genes ex- 
pression in rice leaves showed as genotype- dependent 
fashion [11]. In present work, the specific genes involved 
in plant height and circadian rhythms are matched to 
phenotype bias as well, suggesting that epigenetic regula- 
tion is highly selective and complicated. Limited data ob- 
tained in the present work may have led to finding the 
H3K36me3 modification only on the partial genes that 
were associated with the GCA effect, and these modifica- 
tions did not match the expression level of the genes ana- 
lyzed. The results of H3K27me3 modification profile 
biased to TQ in GL x TQ, but not bias in another Fl hy- 
brids GL X 93-11 suggested that epigenetic modifications 
could be dynamic changes depending on genetic back- 
grounds, developmental stage, environment etc. There- 
fore, more research into epigenetic modifications, such as 
other histone modifications, acetylation and DNA methy- 
lation, are needed to fully understand the genome-wide 
regulation network. 

Extensive studies have revealed that a large number of 
genes exhibited differential expression level between Fl 
hybrids and their parents [8,10,12,13,55,56]. However, 
the regulation of gene expression profiles in the hetero- 
zygous state remains unknown. Our data revealed that 
the positive GCA effect parent in the Fl hybrids could 
determine the phenotypes of the important agronomic 
traits, such as plant height and heading date. Plant 
height and heading date are very important agronomic 
traits that affect rice grain yield [27,57]. Regulation of 
the genes in the signal transduction pathway and synthe- 
ses of metabolic products have been reported to occur 
through SNPs in these gene sequences that could be op- 
timized to the gene expression profiles in elite varieties 
[21]. Our results indicated that the phenotypes of paren- 
tal bias in Fl hybrids could be the consequence of tran- 
scriptome profile or some specific gene expression level. 

One of the most important characteristics of the elite 
parent with high GCA effect is its wide adaptability. The 
wide adaptability of the elite parent could be due to both 
photoperiod and temperature neutral in order to adapt 
different ecological environments. Therefore, rice plant 
could maintain a constant growth period so that it can 
develop higher biomass and grain yield in different eco- 
logical environments. Previous study has reported the 
OsPRR37 had involved in the regulation network of the 
adaptability to grow in different latitudes [45,57]. In present 



work, we found the H3K36me3 modification patterns of 
the OsPRR gene family, OsPRR37 and OsPRR95, showed 
similar patterns in Fl hybrids, high GCA parents 93-11 
and TQ. The results indicated that the specific modifi- 
cation patterns of H3K36me3 in high GCA parents 
might contribute to the adaptability for elite parents 
and their offsprings. It is reasonable to assume that the 
PRR gene family could play an important role in con- 
trolUng the adaptability of an elite parent that is respon- 
sible for GCA effects in rice. 

Conclusions 

The results presented here revealed that the transcrip- 
tome profiles and DEGs in the Fls were highly related 
to phenotype bias to the positive GCA-effect parent. 
The transcriptome bias was also demonstrated by ana- 
lyzing the genes that controlling the specific phenotypes 
of plant height and flowering and H3K36me3 modifica- 
tion on global modification or specific signaling path- 
ways. Our results indicated the molecular basis of GCA 
is both on transcriptome and epigenetic level. 

Methods 

Plant material and phenotype analysis 

Five rice elite varieties, GL, TQ, 93-11, AJ and ZS, were 
chosen for this study. Their pedigree as showed in 
Additional file 6. GCA evaluation was conducted by a dia- 
llele crossing design. Briefly, 10 parallel crosses derived 
from five parents and the phenotypes were collected from 
a field in Wuhan in 2008 and 2009 and from the Hainan 
Island province in 2009. The plots were established in 
triplicate with 30 plants per plot in every season. Three 
agronomic traits, plant height (PH), heading date (HD) 
and grain number per panicle (GNP), were studied. The 
general combining ability (GCA) for each parent was cal- 
culated according to Griffing [58]. 

Nuclear RNA extraction 

The second fully expanded leaves were harvested at the 
secondary branch differentiation stage, immediately frozen 
in liquid nitrogen and stored at -80°C. The leaves from 
the triplicate plots were pooled for RNA extraction. Nuclei 
were isolated from -10 g of frozen leaves using the Plant 
Nuclei Isolation/Extraction Kit (Sigma, St. Louis, MO, 
USA). Total hnRNA was extracted from nuclei using 
Trizol (Invitrogen, Carlsbad, CA, USA) according to 
the manufacturers instructions, and then treated with 
RNase-free DNase I (New England Biolabs, Ipswich, 
MA, USA) to remove any contaminating genomic 
DNA. 

RNA-Seq Library construction 

The niumina mRNA-Seq Sample Prep Kit (lUumina, San 
Diego, CA USA) was used to prepare the sequencing 
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library with 3 \ig of nuclear RNA. Fragmentation buffer 
in the kit was added directly to hnRNA to produce short 
fragments of 200-700 bp, which served as the templates 
for first-strand cDNA synthesis using random hexamers. 
Second-strand cDNA was synthesized followed the 
protocol described in the kit and was purified using a 
QIAquick PGR Extraction Kit (Qiagen, Valencia,CA 
USA) and eluted in elution buffer (EB). The short frag- 
ments were then ligated to sequencing adapters. Suitable 
fragments of approximately 200 bp were selected as 
templates for amplification in a MyCycler PGR instru- 
ment (Bio-Rad, Hercules, CA USA) with the following 
program: denaturation at 98°C for 30 s followed by 15 cy- 
cles of 98°C for 10 s, 65°C for 30 s, and 72°C for 30 s 
plus a terminal hold at 72°C for 5 min. The samples 
were then purified using the QIAquick PGR Purifica- 
tion Kit according to the manufacturers protocol and 
eluted in 30 [iL of EB. One [xL aliquot of the con- 
struct was loaded onto an Agilent Technologies 2100 
Bioanalyzer using the Agilent DNA 1000 Chip Kit 
(Agilent, Santa Clara, CA USA). After verifying the 
size and purity of the DNA fragments, the library was 
sequenced using an lUumina GA II x platform by BGI 
(Shenzhen, China). 



ChlP-Seq library generation 

Chromatin immunoprecipitation (ChIP) was performed 
with antibodies against trimethylated histone H3 on ly- 
sine 27 (H3K27me3, Abeam Cat. #ab6002, Cambridge, 
MA, USA) and trimethylated histone H3 on lysine 36 
(H3K36me3) (Abeam, Cat. #ab9050, Cambridge, MA, 
USA) as described by Saleh et al. [59]. The DNA was ex- 
tracted by adding equal volumes of phenol/chloroform/ 
isoamyl alcohol to each tube and vortexing briefly. The 
DNA was precipitated with 2.5 volumes of 100% EtOH, 
1/10 volume of 3 M sodium acetate (pH 5.2) and 2 [A of 
glycogen (20 mg ml"^) at -80°C. Next, the ChIP DNA 
was used to generate Illumina sequencing libraries fol- 
lowing the manufacturer s protocol of Paired-End DNA 
Sample Prep Kit (Catalog #: PE-102-1001, Illumina, San 
Diego, CA, USA): appropriate fragments approximately 
200 base pairs (bp) in length were selected as templates 
for amplification involving denaturation at 98°C for 30 s; 
15 cycles of 98°C for 10 s, 65°C for 30 s and 72°C for 
30 s; and a final incubation at 72°C for 5 min. The sam- 
ples were then purified using the QIAquick PCR Purifi- 
cation Kit according to the manufacturer s protocol and 
were eluted in 30 \i\ of elution buffer. One microliter of 
the library was loaded on an Agilent Technologies 2100 
Bioanalyzer using the Agilent DNA 1000 chip kit (Agi- 
lent, part #5067-1504). After its size and purity were 
verified, the library was sequenced using the Illumina 
GAIIx platform. 



Analysis of ChlP-Seq reads with the Reads Per Kb per 
Million reads (RPKM) method 

Based on the gene body-specific distributions of both 
H3K27me3 and H3K36me3 [40,60] [10], we normalized 
the ChlP-Seq read counts by calculating the number of 
RPKMs in the gene body region. 

Statistical analysis 

A threshold of more than a 2-fold difference in gene ex- 
pression levels and a False Discovery Rate (FDR) of less 
than 0.05 were used to identify the DEGs. An FDR less 
than 0.001 were used to identify the genes differentially 
modifying histones. The evaluation of the P value and 
FDR of DEGs was performed as described by Audic 
et al. [61] and Benjamini et al. [62], respectively. The t- 
test and correlation analysis were conducted using 
Microsoft Office Excel 2010. 

Quantitation real time PCR 

The confirmation of gene expression level by qRT-PGR 
were described by Song et al. [11], and the methods were 
described by Wang et al. [63] . 

Additional files 



Additional file 1: Figure 51. Phenotypes related to grain yields from 
10 diallele crosses. 

Additional file 2: Figure 52. Plant morphologies of 10 diallele crosses. 

Additional file 3: Figure 53. Plant phenotypes of three parents and 
their Fl hybrids used for transcriptome sequencing. 

Additional file 4: Table 51. The expression levels of the known genes 
involved in the flowering pathway in rice. 

Additional file 5: Table 52. The expression level of the known genes 
involved in GA metabolism and signaling pathways in rice. 

Additional file 6: Figure 54. The pedigree of five parents used in this 
study. The pedigree of five parents was constructed by 153 
polymorphism SSR markers using the Ntsys2.1 software. 
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